LAMMPS (3 Nov 2022)
boundary p p f
kspace_style pppm/electrode 1.0e-7
kspace_modify slab 3.0

include "settings.mod" # styles, computes, groups and fixes
# set boundary in main script because ffield is periodic
units real
# distribute electrode atoms among all processors:
if "$(extract_setting(world_size) % 2) == 0" then "processors * * 2"
processors * * 2

atom_style full
pair_style lj/cut/coul/long 14
# kspace_style and _modify in main script to test different approaches

read_data "data.au-vac"
Reading data file ...
  orthogonal box = (0.285001 -0.368323 -1.043333) to (25.325067 24.671743 43.303354)
  1 by 2 by 2 MPI processor grid
  reading atoms ...
  288 atoms
  reading velocities ...
  288 velocities
Finding 1-2 1-3 1-4 neighbors ...
  special bond factors lj:    0        0        0       
  special bond factors coul:  0        0        0       
     0 = max # of 1-2 neighbors
     0 = max # of 1-3 neighbors
     0 = max # of 1-4 neighbors
     1 = max # of special neighbors
  special bonds CPU = 0.003 seconds
  read_data CPU = 0.007 seconds

group bot molecule 1
144 atoms in group bot
group top molecule 2
144 atoms in group top

# ramping potential difference
variable v equal ramp(0,2)

# get electrode charges
variable q atom q
compute qbot bot reduce sum v_q
compute qtop top reduce sum v_q

# get theoretical charges:
# calculate distance dz between electrodes
compute zbot bot reduce max z
compute ztop top reduce min z
variable dz equal c_ztop-c_zbot

# calculate theoretical capacitance as eps0 * area / dz
variable eps0 equal 55.26349406/10000 # epsilon zero
variable capac equal "v_eps0 * lx * ly / v_dz"

# calculate theoretical charges and deviation of constant potential charges from theory
variable qtheory equal "v_capac * v_v"
variable percdev equal "100 * (c_qtop - v_qtheory) / (v_qtheory + 1e-9)" # avoid divide-by-zero

# constant potential electrodes with ramping potential difference
fix conp bot electrode/conp 0 1.979 couple top v_v symm on
288 atoms in group conp_group

thermo 1
# thermo: step, imposed potential, bottom charge, top charge, theory charge, percent deviation
thermo_style custom step v_v c_qbot c_qtop v_qtheory v_percdev

run 10

CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE

Your simulation uses code contributions which should be cited:

- kspace_style pppm/electrode command:

@article{Ahrens2021,
author = {Ahrens-Iwers, Ludwig J.V. and Mei{\ss}ner, Robert H.},
doi = {10.1063/5.0063381},
title = {{Constant potential simulations on a mesh}},
journal = {Journal of Chemical Physics},
year = {2021}
volume = {155},
pages = {104104},
}
- fix electrode command:

@article{Ahrens2022
author = {Ahrens-Iwers, Ludwig J.V. and Janssen, Mahijs and Tee, Shern R. and Mei{\ss}ner, Robert H.},
doi = {10.1063/5.0099239},
title = {{ELECTRODE: An electrochemistry package for LAMMPS}},
journal = {The Journal of Chemical Physics},
year = {2022}
volume = {157},
pages = {084801},
}
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE

WARNING: No fixes with time integration, atoms won't move (src/verlet.cpp:60)
PPPM/electrode initialization ...
WARNING: System is not charge neutral, net charge = 0.000219 (src/kspace.cpp:327)
  using 12-bit tables for long-range coulomb (src/kspace.cpp:342)
  G vector (1/distance) = 0.16696028
  grid = 8 8 25
  stencil order = 5
  estimated absolute RMS force accuracy = 2.7967012e-05
  estimated relative force accuracy = 8.4221827e-08
  using double precision MKL FFT
  3d grid and FFT values/proc = 2691 672
Generated 1 of 1 mixed pair_coeff terms from geometric mixing rule
Neighbor list info ...
  update: every = 1 steps, delay = 0 steps, check = yes
  max neighbors/atom: 2000, page size: 100000
  master list distance cutoff = 16
  ghost atom cutoff = 16
  binsize = 8, bins = 4 4 6
  2 neighbor lists, perpetual/occasional/extra = 2 0 0
  (1) pair lj/cut/coul/long, perpetual
      attributes: half, newton on
      pair build: half/bin/newton
      stencil: half/bin/3d
      bin: standard
  (2) fix electrode/conp, perpetual, copy from (1)
      attributes: half, newton on
      pair build: copy
      stencil: none
      bin: none
WARNING: Using kspace solver on system with no charge (src/kspace.cpp:314)
Per MPI rank memory allocation (min/avg/max) = 7.725 | 7.774 | 7.823 Mbytes
   Step          v_v           c_qbot         c_qtop       v_qtheory      v_percdev   
         0   0              0              0              0              0            
         1   0.2           -0.030912769    0.030912769    0.031376826   -1.4789807    
         2   0.4           -0.061825538    0.061825538    0.062753652   -1.4789807    
         3   0.6           -0.092738307    0.092738307    0.094130479   -1.4789807    
         4   0.8           -0.12365108     0.12365108     0.1255073     -1.4789807    
         5   1             -0.15456384     0.15456384     0.15688413    -1.4789807    
         6   1.2           -0.18547661     0.18547661     0.18826096    -1.4789807    
         7   1.4           -0.21638938     0.21638938     0.21963778    -1.4789807    
         8   1.6           -0.24730215     0.24730215     0.25101461    -1.4789807    
         9   1.8           -0.27821492     0.27821492     0.28239144    -1.4789807    
        10   2             -0.30912769     0.30912769     0.31376826    -1.4789807    
Loop time of 0.0162257 on 4 procs for 10 steps with 288 atoms

Performance: 53.249 ns/day, 0.451 hours/ns, 616.306 timesteps/s, 177.496 katom-step/s
50.4% CPU use with 4 MPI tasks x no OpenMP threads

MPI task timing breakdown:
Section |  min time  |  avg time  |  max time  |%varavg| %total
---------------------------------------------------------------
Pair    | 0.0010955  | 0.0013022  | 0.0014351  |   0.3 |  8.03
Bond    | 6.13e-07   | 7.795e-07  | 1.031e-06  |   0.0 |  0.00
Kspace  | 0.0034606  | 0.003541   | 0.0036869  |   0.2 | 21.82
Neigh   | 0          | 0          | 0          |   0.0 |  0.00
Comm    | 0.0015915  | 0.0017672  | 0.0019298  |   0.3 | 10.89
Output  | 0.00071528 | 0.00095256 | 0.001506   |   0.0 |  5.87
Modify  | 0.0080243  | 0.0080738  | 0.0081204  |   0.0 | 49.76
Other   |            | 0.0005881  |            |       |  3.62

Nlocal:             72 ave          72 max          72 min
Histogram: 4 0 0 0 0 0 0 0 0 0
Nghost:           1062 ave        1062 max        1062 min
Histogram: 4 0 0 0 0 0 0 0 0 0
Neighs:           6624 ave        6936 max        6312 min
Histogram: 2 0 0 0 0 0 0 0 0 2

Total # of neighbors = 26496
Ave neighs/atom = 92
Ave special neighs/atom = 0
Neighbor list builds = 0
Dangerous builds = 0

Total wall time: 0:00:00
